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ABSTRACT 


Feasibility of utilizing the Random Decrement method in conjunc- 
tion with a Signature analysis procedure to determine the dynamic 
characteristics of an aeroelastic system for the purpose of on-line 
prediction of potential on-set of flutter has been examined . 

Digital computer programs were developed to simulate sampled 
response signals of a two-mode aeroelastic system. Simulated response 
data were used to test the Random Decrement method. A special curve- 
fit approach was developed for analyzing the resulting Signatures. A 
number of numerical "experiments" were conducted on the combined pro- 
cessess. The method was found to be capable of determining frequency 
and damping values accurately from Randomdec Signatures of carefully 
selected lengths. 
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SUMMARY 


An ensemble averaging method for determining the characteristic • 
response function of an aeroelastic system from its turbulence-induced 
random vibrations was developed recently by Henry A. Cole, Jr. The 
most significant feature of this (the Random Decrement) method is that 
very little knowledge about the excitation is required to obtain use- 
ful results (Randomdee Signatures) . Provided an automated numerical 
procedure can be developed to analyze Signatures, the Random Decrement 
method can be very useful for on-line prediction of the on-set of flut- 
ter during subcritical wind tunnel or flight testing. 

The Signature analysis procedure selected for investigation under 
this project is that of a curve- fitting nature. A Randomdee Signature 
is approximated by the theoretical homogeneous solution of the mathe- 
matical model of the aeroelastic system under consideration. An error 
function between the Signature and the analytical expression is defined, 
determined and numerically minimized. Those coefficients in the theo- 
retical solution which lead to a minimum error are said to be the best 
approximations of the dynamic properties of the aeroelastic system. 

Digital simulation techniques were employed to carry out numeri- 
cal experiments for investigating effects of variation of system, data 
acquisition, Randomdee and curve-fit parameters on the accuracy of the 
overall approach. 
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Section 1 


INTRODUCTION 

The need for an accurate, reliable and rapid means to forecast 
the onset of flutter during subcritical wind-tunnel and flight tests 
of aeroelastic structures is widely recognized. That this need exists . 
today is evidenced by the number and variety of current and recently 
completed research and development programs in the United States and 
abroad (References 1, 2 and 3), as well as by the differences among 
methods used by various organizations (References 4, 5, 6 and 7). 

Coupry (Reference 4) and Cole (References 5, 6, 8 and 9) both advocate 
the utilization and on-line processing of wind induced dynamic response 
signals from the test specimen for flutter prediction and failure de- 
tection. Their approaches have obvious advantages over other methods 
in which "controlled" excitations (impulsive, periodic and/or other 
kinds of forces) must be applied to induce responses in those vibra- 
tional modes which eventually flutter. 

Coupry relies on rapid PSD analyses of the response signal and 
modal identification in the frequency domain. The accuracy of the 
method suffers when flutter modes are closely spaced on the frequency 
axis. The method also depends on the PSD of the excitation (in this 
case, the effects of turbulence) being reasonably flat over the fre- 
quency range in which the flutter modes reside. 

Cole first introduced the concept of what was later named the 
Random Decrement method in Reference 5 while investigating applications 
of the correlation functions. He subsequently used it at NASA/Ames 
Research Center for detecting fatigue failures in a Space Shuttle 
wing flutter model (Reference 6). The method received full treatment 
by the inventor in References 8 and 9. 



The Random Decrement method is essentially an ensemble averaging 
procedure which determines the characteristic response function (the 
Randomdec Signature) of an aeroelastic specimen under test from its 
turbulence- induced random vibrations. 

To obtain a Randomdec Signature, one simply collects a number of 
segments of time series representing the random responses of a system, 
and ensemble-averages them. If the system is linear and the excitation 
random, the average time series converges towards the transient response 
of the system due to a set of initial conditions. The order of the sys- 
tem is arbitrary. The initial conditions to which the Signature cor- 
responds can be manipulated almost at will by judicial "triggering" 
(selection of the starting point) of each ensemble member. For failure 
detection and for property identification of nonlinear structures. Cole 
favors triggering at a constant response level. The resulting Signature 
is an approximating time series of the characteristic response function* 
of the specimen in its natural environment at an amplitude which is e- 
qual to the trigger level. For flutter prediction in linear systems, con- 
stant-level triggering is not necessary, and various other methods be- 
come optionally available. Since the ensemble averaging procedure can- 
not be carried out indefinitely in practical situations, the Randomdec 
Signature will contain a certain amount of error which generally makes 
direct (visual) interpretation and determination of system dynamic 
characteristics difficult. The main objective of this study is 
to investigate the feasibility of utilizing the Random Decrement method 
in conjunction with a Signature analysis procedure to determine the 
dynamic characteristics of an aeroelastic system under test for the 
purpose of on-line prediction of potential on-set of flutter. 

*The existence of a characteristic response function for a nonlinear system 
is an assumption. 
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The Signature analysis procedure selected for investigation un- 
der this project is that of a curve-fitting nature. A Randomdec Sig- 
nature is approximated by the theoretical homogeneous solution for the 
mathematical model of the aeroelastic system. An error function between 
the Signature and the analytical expression is defined, determined and 
numerically minimized. Those coefficients in the theoretical solution 
which lead to a minimum error are said to be the best approximations 
of the dynamic properties of the aeroelastic system. 

Digital computer programs were developed to simulate sampled re- 
sponse signals of a two-mode aeroelastic system. Simulated response 
data were used to test the Random Decrement method. A special curve- 
fit approach was developed for analyzing the resulting Signatures. A 
number of numerical "experiments" were conducted on the combined pro- 
cesses. Results of these experiments indicate definite feasibility 
of combination of the approaches. 

Analyses of the Random Decrement method and the curve-fit pro- 
cedures are presented in Sections 3 and 4 of this report, respectively. 
Section 5 deals with the simulation of response signals and the imple- 
mentation of the Random Decrement procedures. Numerical results of 
the investigation are summarized in Section 6. 

Computer programs developed for this study are described in Appen- 
dix A. A sample test case with typical input and output is included 
in Appendix B. 


3 



Section 2 


SYMBOLS 


All symbols used in the text portion of this report are defined 
when they are first introduced. The following is a cross-reference 
of these symbols arranged in alphabetical order. 


SYMBOL DEFINITION 


a. 

i 

A i (s) 

b. 

i 

B. 

l 


B. 

— i 


BC. 


IK 


C. 

1 


l 



0 



ni 


D. 

l 


D. 

— l 


Quasi-steady aerodynamic generalized force in 
the ith mode 

History dependent component of generalized force 
in the ith mode due to aerodynamics 

See definition of C-T^Au) 

Amplitude in the ith mode of the anti-symmetric 
part of the periodic factor of a Randomdec Sig- 
nature 

Best approximation of B. determined by least- 
squares curve-fit 

"Box-car’ T components of the generalized force 
in the ith mode 

Total modal damping coefficient in the ith mode 

Modal damping coefficient of the structure in 
still air 

Modal damping coefficient due to aerodynamics 


Intermediate variable used in the curve-fit pro- 
cess 

Correction factors calculated by the curve-fit 
procedure for the parameters a during the £th 
iteration m 

Amplitude in the ith mode of the symmetric part 
of the periodic factor of a Randomdec Signature 

Best approximation of D. determined by least- 
squares curve-fit 
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e(w) 


E 


E 


l 


Infinite Fourier Transform of measured velocity 
response of an aeroelastic system, including 
effects of low-pass filters 

Error function used to gauge the degree of 
success of a curve-fitting process 

Error function E after & iterations of the curve- 
fit process 


f 

o 


Bandwidth of the response simulation process. 
Hertz 


fjCt) 



g 


Inverse Laplace Transform of the factor in the 
transfer function on the ith mode which distin- 
guishes an aeroelastic system from a purely 
mechanical system. For a purely mechanical 
system, f^(t) is a Dirac Delta function 

Sequences of random numbers used to construct 
in the simulation process 

A dummy integer subscript used in intermediate 
steps in the curve-fit process 


h.(t) 


Impulse response function of the ith mode due 
to an impulse applied at t=0 


hlOO 

h'.’ (t) 
1 


i 


Convolution of f . (t) and hV (t) 
i 1 

Inverse Laplace Transform of the factor in the 
transfer function in the ith mode, which is 
responsible for flutter of an aeroelastic sys- 
tem when its poles move on to the imaginary 
axis 

A subscript identifying modal parameters and 
variables, i=l, .... I 


I 


See definition of i. 


j -Index used tc identify tho-se modes in ..which 

generalized forces are statistically related, 
j=l, . • ■ . J<I 

J See definition of j . 


k Sample counter, k=l,2, . . . K 

k^ Total generalized stiffness in the ith mode 

k.' Generalized stiffness of the structure in still 

1 air 
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k'.' 

X 

K 

L { } 

in 


m. 

x 


M 


m. 

x 


n 


N 


10 


P- 

r mo 



Generalized stiffness due to aerodynamics 

See defnintion of k 

The Laplace Transform operator 

An integer subscript used to associate inter- 
mediate variables in the curve-fit process with 

the various a 's 
m 

Total generalized mass in the ith mode 

Generalized mass of the structure in still air 

Generalized mass due to aerodynamics 

An integer counter used in script form for en- 
semble members in the Random Decrement process. 

It is also used in identifying the time associated 
with each ensemble observation. 

See definition of n 

Average initial velocity in the ith mode after N 
ensemble averages 

Velocity in the ith mode at t=t no , which is also 
the initial modal velocity of the nth member used 
in the Random Decrement process 

Simulated periodic components of the generalized 
forces in the ith mode 


q i (t) 
Co) 
Co) 


Response (displacement) of the ith mode 

Initial displacement in the ith mode as a function of time 
Initial velocity in the ith mode as a function of time 



Mean-squared value of response in the ith mode as a 
function of time 


qfct) 


Expected squared value of response in the ith mode 
as a function of time 


q. Average initial displacement in the ith mode after 

10 N ensemble averages 

q. (t ) Response in the ith mode collected by the Random 

in 11 Decrement process as the nth ensemble member, 

starting at t =0 
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q. 

xno 


Q i (t) 

Q!(t) 



Displacement in the ith mode at t=t no , which 
is also the initial modal displacement of the 
nth member used in the Random Decrement process 

Generalized force in the ith mode 

Part of generalized force in the ith mode 
which is statistically independent of generalized 
forces in all other modes. 

Part of generalized force in the ith mode 
which is statistically related to a similar 
component in the jth mode 

Mean-squared value of the generalized force 
in the ith mode 

Mean-squared value of Q! 

Mean- squared value of QV. 




R. . 
ij 


S N^ 
S (t») 

co x J 


t 


Simulated random components of generalized force 

for the ith mode. Q. is simulated at discrete time points 

ic6 only. 1K 

Intermediate variable used in the curve-fit process 

Intermediate variable used in the curve-fit process 

Coefficients relating generalized forces in dif- 
ferent modes vfhich are statistically dependent 

Laplace Transform variable 

"Signal" paTt of a Randomdec Signature 

Ideal Randomdec Signature obtained by ensemble 
averaging infinitely many statistically inde- 
pendent response samples 

Time 


t’ 


t 

n 


Running variable for time for all members collected 
by the Random Decrement process. It is also used 
as the independent variable for the ensemble average 
of all members (the Randomdec Signature). 

Running independent variable for the nth member 

collected in the Random Decrement process, t =t-t 

n no 
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t 

no 


T 


T 

c 


T. 

1 


(-T./b.) 

x x 


U(t) 


U N (t') 


Time at which the nth ensemble member begins. It 
is also the nth triggering time in the Random De- 
crement process. 

Length of a Randomdec Signature 

Time constants of simulated lovi-pass filters used 
in data acquisition system 

Time constant in the generalized force in the 
ith mode associated with A. (s) 

l 

Real pole of the transfer function in the ith mode 
of an aeroelastic system 

Measurable response of a dynamic system 


Mean- squared value of system response 

Total response function collected as the nth 
member of the ensemble by the Random Decrement 
process 

Ensemble average of (t ' ] 



Value of the kth point in the Randomdec Signature 


y(t') 


2(y k ) 


a. 


a. 


a. 


a. 


a r 


a. 


a_ 


Theoretical homogeneous solution of a purely 
mechanical system used to approximate an aero- 
elastic system 

yCt,) !t'=ks 

A weighting functional 

Equivalent to 
Equivalent to 

Equivalent to w* 

Equivalent to D^ 

Equivalnet to R^ 

I T 

Equivalent to ? 2 W 2 
Equivalent to 
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a. 


Equivalent to D 


2 


a 

m 


General representation of a , a,, . - . . a Q 

L Z o 


Intermediate variable used in the curve-fit process 


2 

Y N 

Y iN (t,) 


6 

Af 

£ 

£ 


"Noise" part of a Randomdec Signature 

Mean- squared value of the noise term in a Randomdec 
Signature 

Expected value of the noise term, Y (t') in the 
Randomdec Signature 1J 

Expected squared value of the noise term in a 
Randomdec Signature attributed to the ith modal 
response 

Sampling period in a data acquisition process 

Frequency resolution of the response simulation 
process in units of Hertz 

Intermediate variable used in the curve-fit process 


Z Larger of £ and 

Z^ Damping factor of the ith mode 

1 

Z^ Real parts of the complex poles of the transfer 

function in the ith mode of an aeroelastic system 


£i 


n 


Best approximation of Z- determined by least- 
squares curve- fit 


K 


£ 

V 


An integer subscript used as a counter for simu- 
lated sampled response data 

An integer superscript used as a counter in the 
iterative curve-fit process 

See definition of p 

V 

Multiplying scale factors for optimize step size 
in searching for least-squares curve fit, V=l,2,3 . 
P v =(o.5)Cv-i] 
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4>. 


<p; 


mgk 

% 

mg 


Integration (dummy) variable for t 

Intermediate variable used in the curve-fit process 

Intermediate variable used in the curve-fit process 


( 0 . 

1 

to! 

i 


Natural frequency of the ith mode 

Imaginary part of the complex poles of the transfer 
function in the ith mode of an aeroelastic system 


to. 


— i 

to. 

1 


Damped natural frequency of the ith mode, to. = 

0 - ??)&i ' 1 

| 2 U l 

Abbreviation for (1 - C. )\j.’ 


to! 

—i 


Best approximation of to. determined by least-squares 
curve- fit 1 
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Section 3 


MATHEMATICAL FOUNDATIONS OF RANDOM DECREMENT 


The mathematical foundation of the Random Decrement method is 
derived in the following for a linear time- invariant system. 


Let t = time, the independent variable, < t < “ ; 

U(t) = measurable response of a system; 

i = a subscript used to identify modal parameters 
and variables, i=l, 2, ... I, (the system is 
of order 21) ; 


q^(t) = modal responses; 


and 


Q^(t) = generalized forces; 

!n(t) = modal response functions after an Impulse applied at t=0. 


The measured (total) response is the sum of all modal responses: 

U(t) = s q (t). 
i=l 


For physically realizable (casual), second-order modes, we have 
h^t) = 0, 0 < t; 

h t (t) = e sin[(l-cf)Vt], t > 0, (1) 

and tu(t) has a finite discontinuity at t=0. In the above expressions, 
is the undamped modal natural frequency and the modal damping fac- 
tor. For any point t (n=l, 2, . . . N) on the time axis, let 

no 

q ino q i Ct )l t =t , 
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mo 


dt 


t=t 


no 


dh (t) 

end h. (t) = , t^O. 

1 dt 


The response function in each mode may he written, in general, 

q.(t) - (P- + 2£.0).q. )h. (t - t ) + q. h. (t - t ) 

H i*- 7 ^nio i 1 - no 7 ino m no 


•A 


COh-ft - x)dx 


no 


where T is an integration (dummy) variable. For each n, we define a 
new independent variable 


t = t - t , 
n no 


and a new modal response function 


Q-- (t ) = q. (t + t ) , 
n h'- n no 75 


which may be written for each n as 


q. (t ) = (p. + 2£.w.q. )h. (t ) + q. _h. (t ) 

An n r ino 1 i n ino i n mo i n 


f tn 

J Q i' 


* I U. (T - t )h(t - T)dT. 
no 7 n 


Since t starts from zero for each n, it is unnecessary to distinguish 

'them m the last "equation. We may, therefore, use a ■conffiori independent 

variable t’ for all t and write 

n 


q. (f) = (p. + 2£.u).q. )h. (f) + q. h. (t') 

An v Ano i lAno i Ano A 


( 2 ) 


/ 


t' 

+ / Q. (x - t )h-(t' - x)dx. 
1 A v no 7 A 7 
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Summing over I modes 3 the response function for each observation be- 


ginning with t ~ t no is 


U (f) = v (p. + 2?.unq. }h. Ct } 

n v J vt ino 1 iino J -* 

x=l 


+ i q. h. (t' ) 
. , hiio l 
1=1 


I /-t’ 


■kl 


+ s 1 Q i ( T - t^h^t' - T)dT. 


The ensemble average of N such observed results is 


1 


N I 


V’ -N S I u n (t') 

n=l i=l 


1 N 1 

= rr 2 2 CP- + 2?.w.q. )h. (f) 

N , . mo l i^ino - 7 ' 
n=l i=l 


i N 1 
+ N 2 , .2 
n=l i=l 


N I p' 

ff * .V Q i tT ' V h i (t ’ - T)dT 

n=l 1=1 n 


1 1 N 

= 2 h . CtM ^ Sp .^] 

1=1 n=l 


i=l n=l 


P « , N 

** h i^4 v- 0 ] 

i=l n=l 


C3) 


I N 


+ N S 2 

i=l n=l* , 0 


/ - t no ,h i tt ’ • T,dT ' 
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1 1 

The quantities —T p. and — 2 q. are the average (over N observations] 

n N , hno N , hno 

n=I n=l 

initial modal velocity and displacement, and respectively. We 

may, therefore, regard the ensemble average of the response functions as 
being of two parts. The first part has a deterministic functional form, 
and its magnitude is established by p^ Q and q shall call this part 

the "signal", and represent it by 

I 

S (t*) = 2 [Cp- + 2£.u).q. )h. (t 1 ) + q. h.(t')]. (4) 

N 1 10 l i n io' i v J “io l J 

The second part, the last term in Equation (3], is the system response 
to the generalized forces, and in our discussions below, it will be 
called the "noise", represented by y^(t']. It is the ensemble average 
of forced responses and its characteristics are dependent on the charac- 
teristics of the generalized forces Q^(t). By interchanging the order 
of summation and integration operations, the noise term in Equation f5] 
maybe re-written in the following form: 

V" 5 = I .2 f h i Ct ' - nr X Qpr - t )]dx 

1=3 o n=l 

Shorthand notations 
N 

Q (t) =2 Q. (x - t ), and 
lN , x i no 

n=l 

Y iN (t ' )s f \Ct' - T]Q iN (T]dr. 

*0 

will be used. 
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For a given set of Q^, the noise term may be evaluated. Since 
Q^(t) may be related to one another, it is necessary to split each 
generalized force into unrelated components Q|(t) and q!J (t) as fol- 
lows: 


Qift) = Q]_00 


J 

+ X R. . Q7 (t) , for each i, 

3=1 13 3 


C5) 


where coefficients R.^ control which and how the generalized forces are 
related . 

We consider first the simple case defined by 

R. . = 0 
i] 

for all combinations of i and j, and where each Q^(t) is an independent 
and stationary random white noise. Then the expected squared value of 
the noise term as a function of t' in each mode may be evaluated by the 
technique outlined in Section 6.2, Reference 10, and is as follows: 


Y- N Ct’) = N Q l f hj(t' - T)dT, 

0 

~2 

where Q. is the mean-squared value of the generalized force Q. (t). 

2 

The expected value of y^pt') is, therefore, 

.T _ _4- t 

-li, / h i (t ' - T)dT - 

1=1 J 
0 

Cross products containing different indices do not appear because of 
the assumed statistical independences among the various generalized 
forces. Each integral on the right-hand side of Equation (6) can be 
evaluated : 


-■15 



f V 2 f V " T ) 2 _ 

J h (f - T)dx - / e sin [C<jJ.(t’ - x)]dr 

n x i 


1 2 2 2 P 

[1 - ?■ - e (1 - £.cos2id.t' + Z-(l - C-) 2 sin2u.t')l 

x 11 1 -i J i 


4£.w. 

i i 


2 V' 

where a>. = (1 - 

l i J i 


For smali values of £ , the above .equation may be replaced by the approximation 


/' 


2 , -2S.w.t’ 

h i Ct ’ - T )dr - (1 - e 1 1 ). 

i l 


(7) 


Using Expression (7) in Equation (6) we obtain 


~2 — , I — j -2£.oi.t' 

7 N Ct') - i 2 QjCl - e 11 )/C4? Ui ). 

i=l 


The last equation may be cleared of the generalized forces explicitly 
by noting that the mean-squared response for each mode is 


q? = lim q'T(t) = Q?/45.w.. 

n i n i x i i l 


Therefore } 


2 1 1 2 
Y\tCT 1 ) * rr 2 qfci - e ) 


N 


N. i 
1=1 


CB) 


and = lim Y*(f) * | \ qj 

t— i=1 
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On account of the assumed statistical independency of the gener- 
alized forces, the modal response components are also independent of 
one another. The mean-squared value of the measured response is, con- 
sequently, the sum of the mean-squared modal responses: 

~2 1 ~2 
U = S q ± . 

i=l 


So that 


Y N = 


io 2 

N 


( 9 ) 


If the generalized forces are not completely independent of one 

another as we have assumed in the first example, the general expression 

of Equation C5) should be used. Instead of Equation (6], a more com- 

2 

plex expression will be obtained for y (t T ). The added complexity, how- 
ever, is only bookkeeping in nature and presents no fundamental dif- 
ficulties. We will demonstrate, via a second example, that the same 
estimate on the noise term is still valid when the generalized forces 
are related to one another. 

Let us restrict our attention to a two-mode system (i.e., I = 2), 
with generalized forces 


Qi - Q" - Q’’, 


and 


Q? - Qo 


In other words, we set 
QjCt) = Q'(t) = o. 


li 




12 


22 
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and R 21 = 0, 


in Equation (5) for 1=2. The noise term in Equation (3) is in this 


Y N (f) 


■i /v 

r\ 


-T)[z qV(t - t no )]dT 


ct ' -x) + h 9 (t' - x)][ .2 q"(x - t )]dx. 


Since Q^(t) and Q^ft) are statistically independent, 

Y N (t,) = k f h i ft ' “ T)dT + W fQ 2 )2 f ^ h i ft ' 

0 0 

= | CCQ”3 2 + W2 )2] f h i (t ' " x;)dT 


+ *«? 2 / 
*|(Q”) 2 / [hjte 


h (t* - x)dT 


T)h 2 Ct' - x)dx. 


x) + h 2 (f - x)] dx 


where (Q") and (Q") are the mean-squared values of Q'.’ (t) and Q” (t) 
1 Z 1 ^ 

respectively. We know from the first example that 


7 7 f t 7 - 7 -2£ja.t' 

[Wp + CQ^ 2 ] / hj_Ct' - x)dx * q j(l - e ) 

J 0 


2 2 _2 ^2 U 2 tI 
h,Ct' - T)dx * q- Cl - e 1 1 ) 


The last integral in Equation (10) is the only one which we have not seen 
until now. Upon carrying out the indicated integration, however, we 
find its magnitude to be of the order £, the greater of and £ . For a 
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lightly damped system, therefore, Equations (8) and (9) can still be 
used to estimate the mean-squared value of the noise term. 

We now summarize results as follows: The ensemble average of N 

sample series of the total response function is the sum of a signal 
term, S^(t'), and a noise term y^(t r ): 

Vt') = S N (t,) + V t,] > 

where the form of S^(t') is deterministic, and its magnitude is controlled 
by the average values of the modal responses and their derivatives at the 
the beginning of all samples (t=t no , n=l, 2, . . . N) ; and where Y^(t’) 
has a mean-squared value which grows from a near- zero value at t'=0, 
to one Nth of the mean-squared value of the measured response at large 
t'. A pictorial representation of U^(t') is shown in Figure 1 for a 
one-mode system. 

If the starting times t of all sample series are selected on 
the basis of the total response, its derivative, or a combination of 
both, in such a manner that S^(t') does not continually diminish with 
increasing N, i.e., if 

lim U (f) = S^Ct’} ^ 0 for all t'. 

N-*» N 

Since S (t 1 ) contains parameters required to specify the system, it 
will be called the Ideal Signature of the system. 

A Randomdec Signature (of N averages) , on the other hand, is the 
truncated function: 

li N (t') = S^ft') + Y N (t'), t r < T, a finite "signature length", 
where S^(t') is defined by Equation (4). 
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The selection of starting times t^ o is called "triggering". 

A proper triggering method is required to guarantee the convergence 
of the ensemble average to a usable* Randomdec Signature in a reason- 
able** time. Two examples of triggering methods are shown in Fi- 
gure 2. With the Level Triggering method, t = every time response 

no 

signal crosses a preselected level, regardless of the sign of the 
slope. With the Zero-Crossing Triggering method, t = every time response 
crosses zero with a plus slope. A third method which is not shown in Figure 
2 is the Every- Point Triggering method where t = every sampling time. 

Each triggering method will lead to a Randomdec Signature with a different 
apparent form. For a linear system, they all contain the same information 
so far as system dynamic properties are concerned. This is because 
-C.w.t' 

h. ft’) = e 11 sinco.t', and 
i i 

. -C-u.t' 

• 1 1 — — — 
h. ft'] = - e ( 5 . 0 ). sina). t' - oucosu.t'] 

i iii ii 

so that U^(t') may be written in the following general form for all 
Randomdec Signatures regardless of the trigger methods used, see Equa- . 
tions 3 , 4 and 5 : 

I -£■(*). t 1 _ 

lh T (t ' ] = 2 e 11 (B-sinoj-t* + D.cosoj.t'] + Y ft'). (12) . 

N ' v ii i i N 

Different triggering methods merely give us different combinations of 
B . and D . . 

i i' 

* The usefulness of a Randomdec Signature depends on the method employed 
to analyze it, and is one of the main objectives of our study. More dis- 
cussions will be found in Section 6. 

**0ur main concern is the minimum amount of data required to obtain a 
useful Randomdec Signature. Time is most appropriately measured in 
terms of the period of the lowest frequency content of interest. 
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Let us turn our attention to the forcing function on an aero- 
elastic structure at subcritical velocities. We will continue to 
assume that the structural properties are invariant with respect to 
time, and that the forcing function is stationary. In other words, 
the following analysis is applicable when the velocity is constant 
and below the flutter velocity. Let us consider the state of affairs in 
a typical mode. The generalized force may be split up into three basic parts 
the first part, Q^(t), is independent of modal responses and is stationary, 
wideband and random; a second part which is dependent on the modal 
response, and can be further divided into aerodynamic inertia, damp- 
ing and stiffness components; and a third part which is dependent on 
the response history. For the ith mode, let 

m! = generalized mass of the structure in still air, 

mV = generalized mass due to aerodynamics, 

m. = m! + mV, 

ill 

c^ = modal damping coefficient of the structure in 
still air, 

cV = modal damping coefficient due to aerodynamics, 

c. = c! + c'.', 

111 

= generalized 'stiffness or the structure in still 
air, 

kV = generalized stiffness due to aerodynamics, 
k. = k! + kV, 

ill 5 

q^(o) = modal initial displacement, and 
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q^(o) = modal initial velocity. 

The equation of motion for each mode of the aeroelastic systems 
can be expressed most readily with the help of the Laplace Transform 
and can be written as 

2 

[nus + c ± s + + A i (s)]q i (s) 

= HipQpCs) + m i q i (o]s + m i q i (o) + c^fo) (13) 

where s is the transform variable, and 


q i (s) = L{q i (t)}, ; 


Q.(s) s L{Q.(t)}. 


The term A^(s) is a general representation of the history-dependent 
components of the unsteady aerodynamic forces. The simplest version 
of A^(s) is a first-order lag: 


A i ( - S ' ) T.s + 1 


where is the quasi-steady aerodynamic force and T^ is the associated 
time constant. 

The following solution of the modal response function is ob- 


tained from Equation (13) : 




[m.Q i (s) 


+ m i q i (o)s + m i q i (o) + c i q i (o)](T i s + 1) 
. - 

(m.s +■ c.s + k.)(T. s +1) + a. 
i i x' x x 


(14) 


The denominator can be expanded, factored and put into the form 

m. (s 2 + 2^!aj!s + u)! 2 )(T.s + b.) 
i xi x JV x ’ 

*The same symbol will be used for the variable and its Laplace Transform. 
The independent variables will always be used to distinguish them. How- 
ever, when a symbol is used to represent an initial condition, the ” 0 ’' in 

parentheses indicates that the quantity in question is a constant (the 
value of the associated variable at t=G or t'=0) 
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with b^ 3 ad and C| related to the original set of coefficients by 
the following expressions: 

b k = 1 - 2 (?[(»[ - qa).}T i3 (15) 

= w?(l + a./k.)/b., a nd (16) 

1 1 x 1 1 

2qw! = [2 E J - (U)’ 2 -- co?)T.]/b.. (17) 

xx xx x i i 

The modal response can now be written as 

Q, (s) + q i (o)s + q-^o) + 2t i UJ i q i (o) 

q i (5) ' 2 — |2 

s x 2 c'-cjJ'-S + co.' 

xi x 

Let 


T.s + 1 
x 


( 18 ) 


T.s + b. 

i i 


L _1 {-y ~ = h"(t), (19) 

s + 2£!io!s + oj ! 1 

11 x 

T 2 4- 1 

L - V; + ~ b > = qct), and (20) 

i i x 


h'.’(t) * f. (t) = h! (t). 

l x i 


(21) 


The time-domain solution for the stationary subcritical response func- 
tion can then be written in terms of h! (t) , the initial conditions, and 
the motion-independent forcing function as follows: 


q^t) = q i (o)h!^(t) + [^(o) + 2^ i u) i q i (o)]hqt) (22) 

+ J Z Q i (r)h!(t - r)dx. 

0 
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As (a^/kp and (T^ad) approach zero. 


b. 

1 


+ 1, 


co! 

i 


-> w. , 

i 


and 2^!o)! 2£.u). , 

11 xi 


so that hV(t) -> h. ft) and f . ft) becomes a Dirac delta. The solution 
i i i 

consequently approaches that for a purely mechanical system. As 
(a^/k^) and (Tuj.. ) deviate from zero, the impulse response function 
h|(t) must be evaluated by the convolution between hV(t) and f^(t). 
The corresponding Randomdec Signature is 


V *’ 5 =. £ Cfeio 

1=1 


+ 2?.o).q. ] h! (t ' ) + q. h!(t’) 

1 l n lO J 1 10 1 


( 23 ) 


+ v N (t-). 

We see that for an aeroelastic system at a subcritical velocity, the 
Random Decrement procedure will produce a Signature which differs 
from that of a purely mechanical system, only in the difference be- 
tween h^(t) and lu (t). In representing the solution by Equations 
(IS) -(22), the above mathematical differences between solutions of 
the two types of systems are more easily related to the physical dif- 
ferences . 

In any case, we see that the addition of the lag term does not 
increase the number of independent initial conditions beyond two 
per mode. 
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Finally, we remark that by (1) setting 1=2, a^aj, T^T^, and 
QjCt) “ Q^ft), (2] linearly combining q^(t) and q 2 (t), and (3) dis- 
regarding all initial conditions, we will be able to simulate with 
Equation (14) the special case used by Houbolt (Equation (49), Refer- 
ence 11). 
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Section 4 


SYSTEM IDENTIFICATION VIA SIGNATURE CURVE FITTING 

The Randomdec Signature is a truncated approximation of the 
characteristic response function (transient response function due 
to initial conditions) of the system in question. It contains in- 
formation on the system characteristics. From subcritical flutter 
testing point of view, it is important to be able to determine the 
natural frequencies and damping ratios from the Signature. The 
method selected to accomplish this objective is described below. 

Our efforts are concentrated on a two-mode aeroelastic system 
which is characterized by the ten parameters wj, w),, Zy Zy , 
a^, a^, and T^ . The generalized masses and m^ cannot be 
determined from Randomdec Signatures because they appear explicitly 
only in forced response solutions. Fortunately their determination 
is not needed in subcritical testing as they are only multiplying 
constants for the responses and do not directly relate to stability. 

To completely specify a Randomdec Signature for a two-mode system, 
four additional parameters corresponding to the initial amplitudes 
and velocities in the Signature as governed by the triggering me- 
thod are required. The data acquisition equipment and process is 
assumed to introduce yet another three parameters: the sampling period, 

6 and two filter time constants, T, and T . The objective of 
the curve-fit procedure is to determine the best estimates for the 
four most important system parameters wj,, and Gi, in the 
presence of all others, from an imperfect Signature in which not 
all effects of the forcing function have been averaged out. We 
will try to curve-fit the Signature, U^Ct'), by 
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I -£lw!t' _ _ 

y(t’) = 2 e 1 1 (B.sinant' + D^osant') ( - 23 - ) 

with u>! = (1 - to! 2 ) 2 to! 

l v l ' i 

In approaching the problem this way, we are assuming that 

(a) In a real test situation, the response transducer and sig- 

nal conditioning equipment and process are reasonably good, so that 
the cut-off frequencies and oj c corresponding to and , respec- 

tively. are very high in comparison with and col,. 

(b) The effects of a^ and are primarily on the natural fre- 
quencies and damping ratios, and by using and in hl(t), most 
of the differences between the impulse response functions of the 
aeroelastic system and a purely mechanical system has already been 
accounted for. Putting it in another way, the convolution indicated 
by Equation (21) is assumed to produce only small differences between 
h"(t) and h^(t). (The consequent error can be estimated by the dif- 
ference between f^(t) and the Dirac Delta.) 

In view of Equation (IS), our approach is valid because flutter 
is reached when or when ^ or £.,) become non-positive. 


We begin the curve-fit procedure by defining an error -measuring 

function E. The most commonly used function for this purpose is the mean- 
squared error: r T 

l I [y(t') - u rt')] 2 dt'. 


In the case where U^(t ' ) is defined only at t'=k6, k=l,2, . . . . K, 
(i.e., if sampled data are used), E is written as 


E ~ K k?i^ k ~ 5 
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where y fc ~ 

and n „ - U. . (k 6) . 
Nk N 


A more general approach is to introduce weighting and use 


E = 


is iz(y k } - zfU^)] 
k=l 


2 


(24a) 


Parameter values B. , D., to! and C! which minimize E will be 
—x —x —x —x 

considered as the optimum approximations of the true parameters.'* 

The method of minimization is based on the Gauss-Newton method (Refe- 
rence 12} modified by Roman (Reference 13} and further modified dur- 
ing the course of this project. An initial estimate of parameter va- 
lues is made. This estimate is upgraded by information based on the 
slopes of the error function E with respect to variation of each of 
the parameters. For purpose of the. following discussion , it will be 
convenient to rename the parameters in accordance with the following: 


a, = 


a. 


44 


cs_ 

3 


= OJ' 


a„ = D„ 


Ur- = B, 


a, - 


44 


a_ 

/ 


- to; 


and a„ = 


V 


(25) 


*The parameters Eh and are not system properties but have to be deter- 
mined in the curve-fit process. 
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Let a°, m= 1,2, 3, 4, 5, 6, 7 and 8 be an initai set of estimated values 

of a , and let 

m 

£ . , . , . £ 

C k = (dz/dy) k 

be the indicated derivatives evaluated at t'-kS during the £th itera- 
tion. In addition, we calculate 


Hi o 

r, = C(3y/3a )5, 
mk k J m k 


0 0 £ 

W = Wgk- 

mg k=1 mgk 

V = e g kf zty d - zC V ! 


K 


0 “ 0 

R 

g k< gk 


2 2 - 12 . 

[<n =- [f .3 [r ), 

L L nng 1 


and 


£ ,£ . £ 
£ = d /a . 
m mm 


In the Gauss-Newton method, if 


]e^ > e = a predetermined band for all m. 


til 

improved, (£ + 1) parameters are obtained by 


$■ * a 1 * db 

m mm 


In our method, an additional optimum step size factor p^, is 
used to improve the rate of convergence and the last expression is 
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modified as follows: 


a<* + V 
m 


£ x A l n 
a + d p . 
m m v 


(26) 


Initially, v = 1, and p^ = 1 is used. 
If the resulting error function 

E C£ + 


the step size is halfed, i.e.. 


P 2 - 0.5, 


and the error function is checked again. The process continues (v = 3, 
4, . . . ) until 

E C * + 1} < e £ . 

The resulting set of new estimates ^ will then be used to repeat 

m 

the entire process until the condition 


a 


(£ + 1 ) 


m 



< £ 


(27) 


is satisfied for all m. 

The curve-fit procedure is outlined in Figure 3. In addition to 
the above. logic, the computer program is implemented to conduct both 
a four-parameteT (one-mode) and an eight -parameter (two-mode) fit each 
time. For the one-mode case, the initial estimate for (a) the natural 
frequency is obtained by counting the number of samples per each '’period” 
of the Signature; (b) the amplitude is the magnitude of maximum peak; 

(c) the damping is set to zero. For the two mode case, the above set of 
initial estimates are used for one of the modes, and the frequency and 
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amplitude results of the one-mode curve-fit are used as initial estimates 
for the other mode. Damping is again set to zero. 
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Section 5 


SUBCRITICAL AEROELASTIC RESPONSE SIMULATION 


The simulation of desired response signals of an aeroelastic 
system may be achieved in the time domain by direct numerical 
integration of the governing differential equations, or by convolving, again 
in the time domain, the forcing function with the impulse response 
function of the system- The simulation of desired response signals may 
also be achieved by complex multiplication of the forcing function* 
with the frequency response - function of the system in the frequency do- 
main. The last approach was chosen in this study because it leads to 
the implementation of a more general computer program which can be easily 
altered to simulate different types of systems. 

The velocity response of an aeroelastic structure in a subcriti- 
cal condition will be simulated. Velocity signal is chosen for 
simulation because in practice it offers an optimum balance between 
love- and high-frequency signal magnitudes under normal test environments. 
From Equation (14) the modal response velocity is 


Pj^s) = sq i (s) - q i (o) 


[ s Q i (s) 


+ sq t (o) - W?q i (o)](T i s + 1) - a-q-^oj/nu 
2 2 

fs + 2Z-0J.S + u>.)(T.s + 1") + a./m. 

^ii i J i ii 


(23) 


The Fourier transform for the modal velocity is obtained when s is re- 
placed by qoj(n -J -l", -oo < u < oo) i n Equation (281: 

’'Including initial conditions which may be considered as special forcing 
functions . 
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2 ® 2 2 
[nwQ C(i))/co + noiq . (o)/ui - q, (o)](nwT + 1) - a q (o)/(i w ) 

z -L JL _L -L JL -L -L -L 

p - (oj) = 

1 2 2 2 
(1 - w /a). + 2qC .(jo/w. ) (pwT. + 1) - a./(m.m.) 

l l i i ill 

The measured total velocity signal has a Fourier transform 

i I 

e (w) = — — 2 Pi M • 

Cn(jaT b + i)CnwT c + i) i=i 

Of course, when e(w) is inverse transformed, the resulting time series 
simulates the velocity response signal. 

The computer program developed during this investigation (Appendices 
A ( 8] uses an existing fast Fourier transform (FFT) subroutine to carry 
out all required forward and inverse transformations. The following re- 
strictions are imposed on the simulation by the FFT*: 

(a) Only equally spaced, sampled response data, in segments 
of finite lengths c a n be simulated, since the FFT is a finite 
discrete Fourier transform procedure. 

(b) The transform of a sampled time series is computed for 
a finite number of equally spaced frequencies only. Both 
bandwidth and frequency resolution are limited. 

(c) The number of data points in the time domain is the 
same as the number of complex frequency components computed 
in the frequency domain. 

Based on procedural requirements of the Random Decrement method, on 
desired frequency resolution and bandwidth, and on computer core size 
limitations, simulated response signals are generated in 2, 048 -point 

^Implications of, and solutions to overcome, these restrictions are dis- 
cussed subsequently. 


( 25 ) 


( 26 ) 
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segments. Displacement and velocity components for each mode at the 
end of a segment are computed and used as initial conditions for the 
following segment- 

The simulation procedure is outlined in the simplified flow dia- 


gram* of Figure 4. Generalized forces are created in the sampled 

time domain first. For each mode, the force can contain a random 

component and two deterministic components. Random components for 

both modes are generated simultaneously from two independent sequences 

of Gaussianly distributed random numbers* F, , in accordance with the 

j K 

following mixing formula: 


Q. = Q. (k6) = z R. .F. 
ik 1J jk 


where i=l,2 and k= 1,2, . . . 2,048 for each segment of simulation, 6 


is the sampling period, and tc <5 are the times at which the simulated sampling 
takes place. Values of the mixing coefficients, R— * axe selected by 
the desired amount of "correlation" between the two generalized forces. 
Sequences and F_,^ have zero means and unity standard deviations. 

Deterministic generalized force components include a periodic 
part, and a "box-car" part BC^.. The amplitudes, frequencies and 

phase displacements of P. , as well as the rise- and fall-times of BC. 
are user selected via input data. The total generalized force for each 


mode is 


Q. = Q. + p. + BC. , i=l,2 . 

IX IK IK IK 

For each mode the above array is transformed into the discrete fre- 
quency domain by the FFT. The transform operates on an array of 4,096" 


*Detailed flow diagrams of all computer programs are found in Appendix A. 
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real points. The first half of the array contains Q. ; the second array 
are filled with zeroes*. After the transformation, only the first 2,048 
complex points in the frequency domain are used. They are combined with 
terms containing initial conditions in accordance with Equation (29), and 
complex-multiplied with modal and filter frequency response functions. 

The real part of the product array of 2,048 complex points is extended 
symmetrically about the frequency origin, the imaginary part antisym- 
metrically. The resulting array of 4,096 complex points is then inverse 
transformed into the time domain. There will be 4,096 real points, of 
which only the first 2,048 are correct data simulating the sampled, fil- 
tered modal velocity response. ; 

The process is repeated for the other mode. Modal responses are 
then summed point by point to obtain the desired total response signal. 

All three series are kept for subsequent Random Decrement processing. 

The following are the two basic relationships among sampling period, 

6 (seconds), frequency resolution, Af .-(Hertz) , and bandwidth, f (Hertz): 

6 (Af) = (1/4096) • 

f = 4096 (A f) 

In addition, at a frequency 2,048(Af), the number of samples per cycle 
is exactly two. 

On account of The tTurricatioTi of the frequency response function, the 
frequency resolutions should be set in such a manner that natural fre- 
quencies of the simulated system is less than f Q /2. -Similarly, 
modal .damping factors causing significant truncation errors in the time 

*A standard preliminary procedure for subsequent convolution via frequency 
domain multiplication, see Reference 14. 
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domain (on the impulse response functions) should be avoided. 

In the computer program, Random Decrement trigger points, t nQ , are 
established by inspecting the total filtered response signals. The 
actual ensemble averaging processes are performed on individual modal 
components as ivell as the total response signal. This represents a 
luxury not enjoyed in real test conditions where modal components 
are not individually available. The added step is used in the simulation 
process to provide study data on how ensemble samples of individual 
responses converge. 
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Section 6 


PARAMETRIC STUDIES 

The digital computer programs were used to conduct a number of 
numerical experiments to establish the feasibility and performance 
characteristics of the Random Decrement, the Signature curve-fit and 
the combined processes. 


6.1 FEASIBILITY OF RANDOM DECREMENT 

The feasibility of applying the Random Decrement method to obtain an 
approximate characteristic response function was demonstrated for a two- 
mode system with and without aerodynamic lag forces.. The two test cases 
were designed to accomplish the goal without relying on the Signature 
curve-fit procedure. 

In the first test case, the aerodynamic lag forces were set to zero; 
the natural frequencies of the two modes are 5,908 Hz and 6.519 Hz (i.e., 
a difference of 9.836%); and the modal damping factors are both 0.020 (of 
critical). Referring to Equations (29) and (30), this test case is specified 
by the following set of parameters: 

Af = 0.02546 Hz 

qi(°) = q 2 (o) = q 2 (o) = 0 

q x (o) = 1.0 in. 

u 1 = 37.12 rad/sec (5.098 Hz) 
u 2 = 40.96 rad/sec (6.519 Hz) 

£ = 0.020 
= 0.020 
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Zero-crossing and every- point trigger methods were used. Results are 
summarized in the series of "log peak plots”* (LPP’s) in Figures 5 
(zero-crossing triggering) and 6 (every-point triggering) . The first 
series of LPP's (Figures 5 -a through 5-f) show the convergence of the 
Signature of the 5.908-Hz Model 3 the second series (Figures 5-g through 
5-1) .the 6.519-Hz mode, and the last series (Figures 5-m through S-q) 
the total system. The number of ensembles used are indicated on the 
figures. These results readily demonstrate that triggering on the total 
response signal does indeed lead to a system Signature which corresponds 
to the sum of two individual Modal Signatures. They also show the left- 
to-right convergence trend characterizing the Random Decrement process, as 
hinted in earlier discussions. The ideal system Signature is shown in Figure 5-r. 

Figures 6-a and 6-b are interim Modal Signatures obtained by the 
every-point trigger method. The amount of response data used to ob- 
tain these interim Signatures are the same as those used to obtain Sig- 
natures shown in Figures 5-d and 5-j, respectively, A comparison of 

*The computer program finds both the positive and negative peaks in a 
Signature; it computes log(Positive Peak Values), and plots the results 
as Peak Nos. 0, 2, 4, ... ; and it computes log(- Negative Peak Values), 
and plots the results on the same graph as Peak Nos. 1, 3, 5, . . Damping 
for a single mode can be determined easily from this plot of the Signature 
peaks by the slope of the curve. For reference, the straight line con- 
nection points (0,0) and (77,-42dB) represents a modal damping of 0.02. 

Figures 5 and 6 also illustrate the print-plot capability developed for 
the Univac 1108. 
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these results indicates that the every-point trigger method requires just 
as much data to. achieve a given Signature accuracy. This method was, 
consequently, deleted from the computer program. 

The following parameter set defines the second test case in which 
modal aerodynamic lag forces are present: 


Af = 0.0200 Hz 

q, Co) = q 2 (o) = q x (o) = q 2 (o) 
= 67.25 rad/sec 
o> 2 = 108.23 rad/sec 


? 1 = 

0.03175 

c 2 

0. 00588 

m l = 

0.0260 

m 2 

0.0520 

R 11 

= 1.920 

R 12 

= 0.385 

R 21 

= 0.960 

R 22 

= 0.000 


2 r 
- /n 

2 , . 


: 2 /i: 

: 2 /i] 

: 2 /i: 
.2 ,, 


a x = 24.69 
a 2 =-93.33 


Ib/in 

Ib/in 
.-3 


T 1 = 3.721 x 10 sec 

T 2 = 3.721 x 10 -3 sec 

T = T =0.000- 
b c 


0 


The case was designed to yield the following aeroelastic system fre- 
quencies and damping 

ojj = 73.54 rad/sec 

t i> 2 = 101.25 rad/sec 

_ 0.00659 

i 

r ' = 0.0338. 

'2 


39 



Simulated responses are shown for this test case in Figure 7. The zero- 
crossing trigger method was used to obtain the set of Signatures in Fi- 
gure 8. Convergence of the Signatures indicated again the feasibility 
of the Random Decrement method.* 

6. 2 FEASIBILITY OF CURVE-FIT 

The feasibility of the curve-fit procedure as a method to identify 
system dynamic characteristics was next established. 

As pointed out earlier in Section 4, the selected study 
approach of this project is to apply a parameter identification techni- 
que which assumes that the Signature is the sum of a free response and 
a small -magnitude forced response of a two-degree-of-freedom system 
without aerodynamic forces. Therefore, checkout of the procedure 
can be logically divided into two steps. The first step is to verify 
that the procedure is accurate if the system to be identified is indeed 
a purely mechanical one. , The second step is to determine whether force- 
fitting the Signature of a system with aerodynamic forces by free 
vibration solutions of a system without such terms would lead to use- 
ful answers for flutter prediction purposes. 

Step 1 was accomplished by adding wideband random noises of vari- 
ous intensities to an ideal Signature to test the capability of the 
curve-fit procedure. The two components of the ideal Signature have 
equal initial amplitudes, and 
0J 2 = 2 05 
C L = 0.018 
- 0.035 

*Detailed analyses of results will be found in Section 6.3. 
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The noise component has a Gaussian distribution, zero mean, and 
RMS magnitudes of up to 18% of the initial amplitude of the Signature. 

One hundred sampled points were used, representing approximately 8 cycles 

of the Signature at the average frequency of the two modes. Near-per- 
fect (within 1%) identification of all properties was accomplished. 

The number of iterations required to determine the least-squares curve- 
fit ranged from 7 to 20. The computer program was consequently modified 
to limit the number of iterations to prevent accidental high-cost com- 
puter runs, in its final form, the program will stop searching for 
the best fit if after 20 iterations, the convergence criterion of Equa- 
tion (27] is still not satisfied. 

It was also determined during this stage of checkout that as few 
as four sample points per cycle of the higher modal response frequency 
are. sufficient for the curve-fit program to identify the correct system 
properties. This feature, together with its ability to determine in- 
dividual modal damping values when the natural frequencies are very 
close* (and, consequently, with the total Signature displaying the char- 
acteristic "beating" of Figures 5-m through 5-f) are among the princi- 
pal advantages of the procedure. In all subsequent investigations a 
sampling rate corresponding to approximately six points per cycle of the 
average frequency was used. 

Figure 9 is a typical Simulated Randomdec Signature presented to 
the curve-fit routine. Figure 10 is the curve which fits the Signature 
in Figure 9 in the least-squares sense. Figure 11 shows the convergence 
paths of a typical parameter and the mean-squared error function. 


*This ability was verified via test cases to be described later. 



6.3 FEASIBILITY OF THE COMBINED PROCEDURE 

Subsequently, the feasibility of combining Random Decrement and curve- 
fit procedures was demonstrated by curve fitting a Randomdec Signature 
of suitable length* obtained for the test case on page 39. Zero-velocity 
properties, i.e., those of the mechanical system, are 
w = 67.25 rad/sec 


u) 2 = 108.23 

rad/ sec 

= 0.03175 


^2 = 0,00588 


m^ = 0.0260 

lb-see 2 /in 

- 0.0520 

2 

lb-sec /in 


_3 

A lag time constant of 3.721 x 10 sec, or exactly one fourth of the 
natural period of the first mode, was assigned. The quasi-static 
aerodynamic forces are 

a^ = 24.69 lbs/in 

a = -93.33 lbs/in 

corresponding to 

2 

a./k. = a, /mm = 0.2100 
11 111 

2 

a 2 /k 2 = = 0.1422 

The selected characteristics for the aerodynamic forces will simulate 
conditions for a velocity close to flutter**. The resulting values of b^ 
and b^ are 1.012 and 0.979, respectively. The factors (T^s + 1 3 / (h s + b^) 
in Equation CIS) will, consequently, not change the form of the impulse 
response function of the system. Based on Equations (15), (16) and (17), 
the calculated natural frequency for the low-frequency mode will increase 

*The effect of Signature length on curve-fit accuracy will be investigated 
in Section 6,4 

**So far as damping values are concerned. The natural frequencies were 
deliberately kept separated by a significant amount to minimize the. cost 
of the computer run. 
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from a zero-velocity value of 63.25 rad/sec to 73.54 rad/sec, while 
that for the other mode will decrease from 108.23 rad/ sec to 101.25 
rad/sec. The corresponding changes in damping are: 

(e = 0.03175) -> (e| = 0.00659) 

C? 2 = .0.00588) •+ (z' 2 = 0.0338) 

The zero-crossing trigger method was used to obtain the Ran- 
domdec Signatures. A Signature length of 50 samples was used. The Signa- 
tures are shown in Figure 8. The following is a summary of numerical 
results from the test run: 



No . of Numerical Re sults of Combined Pro ce ss 

Ensembles; Fr eq £ Error, rad/sec Damping a nd Err or 


( Averaged 

“1 ! 

Error 

“2 

Error 

4 

Error : 
1 

Error 

+ i 

234 

j 

73.48 

+0. 06 

103.91 

+2.66 

.00806 

+.00147 

.0541 +.0203, 

1 461 

I 

73.54' 

0.00 

102.17 

+0.92 , 

.00708 

+.00049 

. 0373 + . 0035; 

914 

73.61 

+0.07 

101.42 

+0.17 

.00719 

+.00060 

.0386! +.0048; 

1806 

73.53 

-0.01 

101.62 

+0.37 ! 

.00590 

-.00069 

.0426' +.0088. 

j 

3612 

73.47 

-0.07 

i 101. 24 

'-0.01 

.00524 

-.00135 

.0351 +.0013! 

Exact 

73.54 

101.25 

.00659 

j. 0338 “| 


The exceptionally good agreement between the exact and numerically obtained 
results indicate that the influence on the aerodynamic forces is, 
indeed, primarily on the frequencies and damping values and that the 
impulse response function of the aeroelastic system can be accurately 
approximated by the use of those for two second-order modes. The results 
also indicate that, when applied properly, the combined process is truly 
a useful tool for flutter prediction with subcritical test data. 


6.4 SIGNATURE LENGTH AND AVERAGING TIME 

* 

Effects of Signature length on accuracy are demonstrated via results 
of two series of tests described below. Problems to be expected in selecting 

"^Obtained with a CDC66Q0 version of the developed computer programs. 
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a suitable Signature length for curve-fitting. for a system whose dynamic 
properties are totally unknown will be discussed in Section 6.5. 

In the first test series, parameters are set up to specify a system 
which may be described as having separated natural frequencies (by +_ 15% 
of the average frequency) and moderate damping (0.04 of critical for each 
mode) . Case parameters are 
Af = 0.020 Hz 

q x (0) = q 2 (0) - q x (0) = q 2 (0) = 0 
a) = 73.91 rad/sec 

(^2 = 99.55 rad/sec 

= 0.040 
C 2 = 0.040 

m = 0.0260 lb-sec 2 / in 

1 

m 2 = 0.0520 lb-sec 2 /in 

R.^ = 1.920 lb/ (lb-sec 2 /in) 

R ^2 = 0.385 lb/ (lb-sec^/in) 

R 21 = 0.962 lb/(lb-sec 2 /in) 

R 2? = 0.000 



Signature lengths of 12, 25, 37, 50, 75 and 100 samples were curve-fitted. 
For reference, the number of samples per cycle at the higher natural fre- 
quency is approximately 5, and the period of beating between two sine waves 
with frequencies equal to and would be approximately 40 sample periods 
Numerical results indicate that the natural frequencies are relatively easy 
to determine except when the Signature length was too short (the 12- or 25- 
point case where the curve-fit program did not have sufficient data to pro- 
cess) , or when it was too long (the 100-point case, due to the excessive 
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weighting of the tail section of the Signature which contains mostly 
random forced response signals}. In all other cases, natural frequencies 
were accurately determined (by within 2.5% of the exact value} after the 
number of ensembles used for averaging reaches 1,000. 

Damping determination, on the other hand, requires considerably 
more data- -approximately 2,000 ensembles were required to determine 
damping to within 0.007 (out of 0.040} , for the 37-, 50-, and 75-point 
cases. The corresponding average time is approximately 145 seconds. 
Damping accuracies are summarized in Figure 11 . Based on these test 
cases, it may be said that Signatures shorter than a beat period are 
usable for the two mode case. On the other hand, a Signature length 
greater than 200% of the beat period requires too much averaging time to 
be considered practical. An optimum Signature length appears to be the 
50-point case which is 125% of the beat period. 

In the second test series, natural frequencies are spaced much closer 
(within +_ 1.2% of the average frequency) and the damping values are 0.04- 
for one mode and 0.005 for the other. The case parameters are: 


u = 84.72 

rad /sec 

u 2 = 86.73 

rad/ sec 

Z 1 = 0.040 


C 2 = 0.005 



all other parameters are the same as the previous test case. 

The beat period of two sine waves at frequencies equal to oi and ui_ 

I ^ 

is 512 sample periods. Signature lengths of 25, 50, 75, 100 and 512 
were curve-fitted. Both the 25- and the 50-point Signatures failed to 
yield useful results after 16,836 averages. For the 75-point Signature, 
the natural frequency prediction did not yield useful results* until 


♦When the natural frequencies are close to each other, the accuracy 
criteria should be based on the measurement of the difference of the 
natural frequencies, rather the the natural frequencies themselves. 
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after 8,418 averages. For the 100-point case, 4,202 averages are required, 
for the 200- and 512-point case, 2,103 averages. The optimum Signature 
length for both natural frequency and damping predictions for this sample 
test series was found to be the 200-point Signature after 2,000 averages, or 
about 150 seconds at the 85 rad/sec average frequency. Damping measurement 
results are plotted against the number of averages for Signature lengths 
of 100, 200 and 512 in Figure 13. 

6.5 APPLICATION NOTES 

The practical usefulness of any on-line flutter prediction process 
greatly depends on the amount of data it has to use at near critical 
velocities. The ideal situation is ttfhen the process requires so little 
data that the need to hold constant velocity at various subcritical 
stages is eliminated. With respect to the overall approach adapted in 
this investigation, this consideration leads to the requirement of deve- 
loping a method to select the proper Signature length So that the Random 
Decrement averaging time can be minimized. 

The following factors must be considered: 

(a) Based on discussions in Section 3, the rms value of the residual 

forced response term (the noise term starts from zero for t’=0, 

and increases with t' in accordance with Equation (8). The expected 
signal -to-noise ratio of the Randomdec Signature for a given number of 
averages is, therefore, maximum for small values of t’. Consequently, 
using a shorter length of the Signature would require less averaging and 
would yield answers with given expected , accuracies in a shorter time. 

(b) For a system known , a priori, to have only one degree of freedom, 
the Signature length and averaging time required to determine damping 

are derived in Reference (9) . 
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(c) For systems with two degrees of freedom, the beating of the 
two Modal Signatures makes it necessary to select a Signature length 
which, in accordance with results described in Section 6.4, is some- 
where between 50% to 125% of the beat period. Once the proper Signature 
length has been selected, the curve-fit procedure can be used with 
confidence . 

This suggests the following approach in applying the techniques 
put together during this investigation. It is assumed that (a) on-line 
Random Decrement and curve-fitting programs (or hardware) have been 
implemented, (b) real-time display of the developing Randomdec Signature 
is available, and (c) it is possible to select the length and sample 
density of the Signature to be presented to the curve-fit program. Starting 
at a constant low velocity where the danger of explosive flutter is not 
present, obtain an accurate reference Randomdec Signature by using a 
large number of averages. A high sampling rate (about 16 points per 
cycle of the response signal) should be used so that the beat period can 
be determined from the visual display. Needless to say, the Signature 
should be long enough to cover more than one beat period. Select from this 
long and dense Signature the proper length (about one beat period) and 
sample density (about 4-6 points per cycle) and present it to the curve- 
fit program. With the properly selected Signature, we are assured of 
fairly accurate frequency and damping results. The minimum number of 
ensembles required to produce frequency and damping data to a given 
accuracy with respect to the reference Signature is then determined 
exp erimentally. 
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The velocity sweep can then begin. The real-time display of the 
current Randomdec Signature and results of the previous curve-fit 
analyses will be used to select the Signature length -used for curve-fit, 
which will be actuated as soon as a sufficient number of ensemble aver- 
ages have been reached. 

It is, of course, possible to automate the above procedure and use 
preliminary curve-fit results to compute the proper Signature length 
for a second, more refined Signature analysis. 
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Section 7 


CONCLUSIONS AND RECOMMENDATIONS 


7.1 CONCLUSIONS 

We have demonstrated analytically in Section 3 that applying the 
Random Decrement data reduction procedure to the random responses of a 
multi-mode system will produce a Randomdec Signature which is an approx- 
imation of the characteristic response function of the system. The 
error contained in the Randomdec Signature is in the nature of a resi- 
dual forced response ivhose mean-squared value decreases as one over the 
number of ensemble averages . 

Using computer programs developed in this project, sampled random 
responses of a two mode system were simulated. The feasibility of the 
Random Decrement procedure was established via a number of numerical 
experiments on different simulated two-degrees -of-freedom systems, including 
both a purely .mechanical and an aeroelastic system. 

A parameter identification procedure using least-squares curve- 
fitting of the Randomdec Signature was adapted. The method was found 
to be capable of determining frequency and damping values accurately 
from Randomdec Signatures of carefully selected lengths. 

For optimum results, a Signature length between 50% to 125% of 
the beat period created by the two frequencies of the Modal Signatures 
should be selected. The number of ensembles required to produce accurate 
damping results by the combined process is found to be approximately 2,000. 

The study was limited to one- and two-mode systems. 
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7.2 RECOMMENDATIONS 


Simulation studies in this project were limited in scope on ac- 
count of computer costs. Since the programs for response simu- 
lation, Random Decrement and curve-fit analysis can all be implemented 
on less powerful computers, exhaustive statistical studies are recom- 
mended using such computers for reasons of economy. 

While flutter usually occurs on account of interaction between two 
modes of an aeroelastic system, the presence of other non-flutter inodes 
at nearby frequencies cannot be denied in practical situations. The 
Randomdec Signature will contain characteristic responses in all modes. 
The Signature analysis (curve-fit) program should, therefore, be extended 
for such situations. On-line flutter prediction programs or hardware 
can then be implemented by suitable modification of the application ap- 
proach suggested in Section 6.5. 
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figure 4 Simplified Block-Diagram Showing Response Simulation Process 
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Figure 5-m System Signature (131 averages, 21 sec) 
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Figure 5*1 System Signature (261 averages, 42 sec) 
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Figure S-o System Signature (540 averages, 84 sec) 
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Figure 5-p System Signature (1062 averages, 168 sec) 
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Figure 5-q System Signature (2108 averages, 336 sec) 
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Figure 5-r Ideal Signature 
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Figure 10 Least-Squares Curve-Fit Result for Signature of Figure 9 
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Figure 12-b Accuracy of Damping Determination vs. No. of Averages 
(Signature Lengths of 12, 25 and 100 samples) 
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SUMMARY 


Program RADCUF synthesizes response data of an aeroelastic 
structure, generates Randomdec Signatures and determines the system 
dynamic characteristics (frequency and damping for each of the two 
modes of the system) via curve-fit procedures. The program is used 
to conduct parametric studies, and to verify the practicality of both 
Random Decrement and curve-fit procedures via numerical experiments. 

The program is written in Fortran IV and is compatible with 
the UNIVAC 1108 Exec 8 version compiler. .The overlay structure of 
this program as implemented in the UNIVAC computer is shown in Fig- 
ure A-l. 

Data-flow during execution can be summarized as follows: 
The main program, DRVR, resides in the primary link and directs the 
logical flow of four subprograms of the secondary link. The first 
subprogram called by DRVR is IDENT which provides the identification 
(name of the program, bin number , run. ID, job number, date and time) 
of the job on output plots. The second subprogram called is FLTR 
which generates random responses of a simulated aeroelastic system 
and obtains Randomdec Signatures from the data. The third subpro- 
gram called is CURVFT which curve- fits the Signature generated by 
FLTR. The last subprogram called is ENDJOB which writes an END OF 
JOB on the last page of the output plots. 



DRVR 


PRIMARY 

LINK 



SECONDARY CONTROL PARAMETERS 

LINK § SUBPROGRAMS CALLED 


r " As - r-Tk 



FIGURE A-l OVERLAY STRUCTURE -PROG RAM RADCUF 
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DESCRIPTION OF PROGRAM 


2-1 Primary Link 

The primary link is named DRVR. The primary function of 
this program is procedural only, and it does not contain executable 
statements. 

2-2 Secondary Link 

Two major programs reside in the secondary link. They are 

FLTR and CURVFT. Their executions are mutually exclusive, and each 
of these programs can be executed independently from the other. Al- 
though program CURVFT normally requires input data generated by pro- 
gram FLTR, this data could be written on Tape Unit 15 by any other 
program, so long as the format is compatible. 

2-2-1 FLTR 

FLTR is a program that resides in the secondary link of 

the system. It creates an environment for the input of general con- 
trol parameters, generates response data of a two-mode aeroelastic 
system and computes Randomdec Signatures. Basic functions of this 
program are; 

(a) Accepts input data. 

(b) Generates response data. 

(c) Calls subroutine FPLT to plot intermediate or final 
results . 

(d) Calls subroutine PKPLT to (1] determine and plot the 
peaks of Randomdec Signature, (2] determine initial 
estimates fdr subsequent curve-fitting and (3) esta- 
blish initial conditions for subsequent response data. 



(e) Calls subprogram RNDMC1 to generate Randomdec 
Signatures . 

(f) Writes final Randomdec Signatures, initial para- 
meter estimates, and all system and Random Decrement 
parameters on Tape Unit 15. 

2-2-2 CURVFT 

CURVFT is a curve-fitting algorithm that operates on Randomdec 
Signatures created by FLTR. The following is the procedure: 

(a) Reads Tape Unit 15 for Randomdec Signatures and 
all other data written by FLTR. 

(b) Sets up an environment to curve-fit the Randomdec 
Signature with a four-parameter expression corres- 
ponding to a one-mode approximation of the Signa- 
ture, 

Cc) Calls subroutine CF which accomplishes the actual 
curve-fitting procedure. 

(d) Calls subroutine FPLT to plot the Randomdec Signa- 
ture being curve-fitted, the analytical expression 
after convergence of the procedure and convergence 
paths of all parameters, and the error function. 

All problem oriented parameters are tabulated by 
FPLT also. 

(e) Repeats steps fh) , Cc], and (d) above for an eight- 
parameter curve- fit. 

ff) Calls subroutine CF to provide the final output of 
the program under simple English text headings. 


A- 4 
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PROGRAM INPUT 


Control parameters input to program RADCUF are made via 
cards utilizing Logical Unit 5 in the UNIVAC 1108 system. Four cards 
are used to define one typical case. The following logic is designed 
in the program for multi-case execution. 

In the first execution, the program will read the first set 
of four data cards and. execute until the generation of the Randomdec 
Signature for this case is completed. Then the flow of execution will 
be directed to read another card corresponding to Card Form 4 or to 
find an EOF (End of File). If a Card Form 4 is encountered, the pro- 
gram will execute utilizing the new parameters given by this card 
and previous input parameters given by Card Forms 1, 2 and 3. If an 
EOF is encountered, the program then will expect to find either another 
EOF or a completely new set of four data cards. If the second EOF is 
encountered, the program then proceeds to curve-fit the generated data. 

The following is the card stream for typical multi-case type 

execution: 

CARD FORM 1 
CARD FORM 2 
CARD FORM 3 
CARD FORM 4 

CARD FORM 4 

CARD FORM 4 

CARD FORM 4 

EOF 

CARD FORM 1 
CARD FORM 2 
CARD FORM 3 
CARD FORM 4 

EOF 

EOF 

A-S 



It should be noted that Randomdec Signatures for all cases to be 
investigated are recorded on Tape Unit 15 before curve-fitting 
begins. 

Input Card Form 1 

Parameters: DF, MU(1), MU(2) , QQ(1) , QQ(2), PCI), P(2), A(l), A(2), 

LD (I) , LW(1), LD (2) , LW(2) 

Format: 9F6.0, 416 

Input Card Form 2 

Parameters: KREF, NREF, NPRT, METHR, NRR, METHS, NRS, NUT, STN, RQ, 

AV (1) , AV (2) , SD (1) , SD (2) 

Format: FS.O, 214, 2(12, 13), 5X, II, 6F8.0 

Input Card Form 5 

Parameters: NA(I), NA(2) , NB, NC, ISRS, IERS, ISC, IEC, ISS, IES , 

ISF, IEF , ISP, IEP , NSF , ISYM, NPF 

Format: 4F7.1, 1314 

Input Card Form 4 

Parameters: NU(1) , NU(2) , Z(l), Z(2), R(l,l), R(l,2), R(2,l), R(2,2) 

C(l) , C (2) , NP (1) , NP (2 ) , TH(1), TH(2) 

Format: 2F4.0, 12F6. 0 

All input parameters are defined in Table A-l. 
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TABLE A-.1 (Part 1) 


TEXT 

SYMBOL 

FORTRAN 

NAME 

DEFINITION 


UNIT . 

RANGE 




6 


sampling period 


sec 

2it/ (4096 (DF) 

Af 


frequency resolution of response simulation process 


Hertz 


2tt CAf ) 

DF 

frequency resolution of response simulation process 


radians/ sec 


m i 

MU (1) 

generalized mass of Mode 1 (including aerodynamic effects) 


mass 


m 2 

MU(2) 

generalized mass of Mode 2 (including aerodynamic effects) 


mass 


qjCo) 

QQU) 

initial displacement of Mode 1 


length 


q 2 (°) 

QQC2) 

initial displacement of Mode 2 


length 



P(i) 

initial displacement of Mode 1 


length/sec 


q 2 (°) 

P(2) 

initial displacement of Mode 2 


length/sec 


a i /C27TA£)" ! 

A Cl) 

magnitude of aerodynamic lag force per unit displacement in 
Mode 1 


mass/ length 


a ? / (2irA£) 2 

A(2) 

magnitude of aerodynamic lag force per unit displacement in 
Mode 2 


mass/length 



LD (1) 

rise-point no. of Box Car part of generalized force in Mode 1 
rise point = (LD (1 ) ) C5 ) seconds after initiation of simulation 

5 

1-2047 


LW(1) 

fall-point no. of Box Car part of generalized force in Mode 

1 

5 

2-2048 


LD (2) 

rise-point no. of Box Car part of generalized force in Mode 

2 

S' 

.1-2047 


LW(2) 

fall -point no. of Box Car part of generalized force in Mode 

2 

6 

2-2048 


KREF 

(Random Decrement trigger level) / (rms response) 


non-dimensional 



NREF 

number of points in each segment of simulated response data 
be used for Random Drecrement processing 

to 

6 

2-1536 


NPRT 

number of points required in Random Signature 


5 

1-512 
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TABLE A- 1 (Part 2) 


TEXT FORTRAN 


SYMBOL 

NAME 

DEFINITION 

UNIT 

RANGE 


METHR 

trigger method [METHR > 0) = level trigger method 

[METHR = 0) = zero-crossing-with-positive 
slope trigger method 

[METHR < 0) = modified zero-crossing trig- 
ger method 

non-dimensional 



NRR 

number of response points to be simulated/2048 

non-dimensional 



METHS 

not used in final version of program 




NRS 

not used in final version of program 

non-dimensional 



STN 

signal-to-noise ratio [maximum response)/ (std . dev. of noise) 

non-dimensional 



RQ 

external key for random number generator 

non-dimensional 



AV (1) 

mean value of random number sequence used to generate random 
component of generalized forces 

non-dimens ional 



AV(2) 

mean value of random number sequence used to generate random 
component of generalized forces 

non-dimensional 



SD(1) 

standard deviation of random number sequence used to generate 
random component of generalized forces 

non-dimensional 



SD (2) 

standard deviation of random number sequence used to generate 
random component of generalized forces 

non-dimensional 


1/ C2ttT 1 (Af)) 

NA ( 1 ) 

cut-off frequency of generalized force in Mode 1 corresponding 
to aerodynamic lag forces 

2tt [A f ) 

1-4096 

i/( 2 rr 2 (Af)) 

NA (2) 

cut-off frequency of generalized force in Mode 2 corresponding 
to aerodynamic lag force 

2ir(Af) 

1-4096 

l/(2TrT b CAf)) 

MB 

first-order low-pass filter cut-off frequency 

2ir[Af) 

1-4096 

l/C2irT c CAf)) 

NC 

first-order low-pass filter cut-off frequency 

2tt [Af ) 

1-4096 


TSRS 

starting point number of transform of forcing function to be 
plotted 

2tt(A£) 

<IDRS 


IERS 

ending point number of transform of forcing function to be plotted 

2tt (Af ) 

<NPRT 



TABLE A- 1 (Part 3) 


TEXT 

SYMBOL 

FORTRAN 

NAME 

DEFINITION 

UNIT 

RANGE 


ISC 

starting time (point no.) of response data to be plotted 

6 

<IEC 


IEC 

ending time (point no.) of response data to be plotted 

<5 

<NPRT 


ISS 

starting time (point no.) of Signatures to be plotted 

6 

<IES 


IES 

ending time (point no.) of Signatures to be plotted 

<$ 

<NPRT 


IEF 

if set >0, interim Signatures will be plotted 

(plot controller) 

0 or positive 


ISF 

not used in final version of the program 




ISP 

not used in final version of the program 




IEP 

not used in final version of the program 




NSF 

number of subframes to be plotted per page 


1,2 or 3 


ISYM 

plot controller, set to zero for print plots, to 35 for 
SC4020 plots 

(plot controller) 

0 or 35 


NPF 

number of horizontal divisions per subframe 

1 + function plotted 


w i /(2u(Af)) 

NU(1) 

natural frequency of Mode 1 

2tt (Af) 

1-2048 

w 2 /(27r£Af)) 

NU (2) 

natural frequency of Mode 2 

2ir (Af) 

1-2048 


Z(l) 

damping factor of Mode 1 (fraction of critical modal 
damping) 

non-dimensional 


? 2 

Z(2) 

damping factor of Mode 2 (fraction of critical modal 
damping) 

non-dimensional 


R u m i 

R(i,D 

magnitude and mixing coefficients of random components 
of generalized forces 

2 

mass -length/ sec 


R 12 ra i 

R(l,2) 

magnitude and mixing coefficients of random components 
of generalized forces 

2 

mass -length/ sec 



•^ 3 > 
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TABLE A-l (Part 4) 


TEXT 

SYMBOL 

FORTRAN 

NAME 

DEFINITION 

UNIX 

RANGE 

* 21*2 

R(2,l) 

magnitude and mixing coefficients of random components 
of generalized forces 

2 

mass-length/sec 


R 22 nl 2 

R(2,2) 

magnitude and mixing coefficients of random components 
of generalized forces 

2 

mass-length/sec 



ccn 

amplitude of sinusoidal component of generalized force 
in Mode 1 

? 

mass-length/sec 



C(2) 

amplitude of sinusoidal component of generalized force 
in Mode 2 

2 

mass-length/sec 



NP ( 1 ) 

frequency of sinusoidal component of generalized force 
in Mode 1 

2 

mass-length/sec 



NP (2) 

frequency of sinusoidal component of generalized force 
in Mode 2 

2 

mass-length/ sec 



TH(1) 

phase displacement of sinusoidal component of generalized 
force in Mode 1 

radians 



TH(2) 

phase displacement of sinusoidal component of generalized 
force in Mode 2 

radians 



NUT 

Units to be written on Signature plots 
(NUT = 1) = in/sec; (NUT =2) = cm/sec; 
(NUT = 3) = ft/sec; (NUT = 4) = m/sec; 
(NUT = 5) = Rad/sec 





4 


FLOW CHARTS 


4-1 FLTR Procedure 

The diagram below shows the order of execution in the response simulation program, 
FLTR. Contents in the' lettered blocks in the- diagram are described in the cor-, 
respondingly lettered descriptive paragraphs following these flow charts. 
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B. 


C. 


B. 


E. 
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ir 



l DETERMINE TERMINAL MODAL 
F | VELOCITIES AND 
' I SAVE AS INITIAL 
j VALUES FOR NEXT SEGMENT 
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INITIALIZE INRS 


l' NORMALIZE THE SUM 
j SIGNATURE BY THE 
TOTAL NO. OF 
TRIGGERS (NTK) 



PKPLT 


MITE FINAL > 
NORMALIZED SIG, 


UNIT 15 



1 


WRITE FINAL v 
SIG ON PRINTER. 


CONTROL PARAMETERS 
AND FINAL NORMALIZED 
SIGNATURE 


RETURN TO READ 
NEXT FORM-4 FROM 
THE INPUT STREAM 
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PROGRAM STOP 25 MAY OCCUR IF DECK 
IS NOT ARRANGED IN CORRECT ORDER 



PROGRAM STOP 35 MAY OCCUR IF DECK 
IS NOT ARRANGED IN CORRECT ORDER 



I RETURN TO BEGINNING OF PROGRAM 
/ TO START AN ENTIRELY NEW CASE 




A-15 





Block A . JRQ is an internal variable, which contains the key for the 
generation of the next random integer by NRAND. It is initialized 
to the user specified input value, RQ, each time a new case is started 
with respect to reading a Card Form 4. 

INRR is an integer program counter which indicates the num- 
ber of segments of response data which have been processed for a given 
case. When it exceeds the required number of segments given by the 
user as NRR, "End-of-Case" procedures are initiated. 

The variable NTM is 0-origin integer sample point counter in 
the time domain. The first point of the second response segment is con 
sidered point No. NTM = 4096 and etc. 

NTK Is an integer to be assigned as the total number of sum- 
mation of input signature. 

AK, BK and CK are constants to be assigned for the generation 
of complex function. 

Fortran arrays SRF (L) , SRG(L), S(L) for L=l, 4096, accumulate 
Signatures (triggered responses) of Mode 1, Mode 2 and the sum of Modes 
1 and 2, respectively. 

The variables P(l) and P(2) are the initial values of veloc- 
ities of Modes 1 and 2. respectively for the first segment of response 
simulation. 

The variables QQ (1 ) and QQ(2) are the initial displacements 
of Modes 1 and 2, respectively, used by the program in generating 
the first segment of the simulated response data. 

The initial velocities and displacements of subsequent seg- 
ments of simulation are obtained according to explanation for Block 
F below. Simulation will be started from initial 
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input values each time a new case is started after reading a Card 
Form A, 

Block B . The arrays (QR1 (L) and QR2(L), L=l, 4096, are used for storage 
of pseudo-random number subsets used for the generation of forcing func- 
tions. 

The Fortran function NRAND is used to generate the next pseudo- 
random integer from the previous one. This integer, converted to float- 
ing form and stored in the first location of the mode array, is used as 
a trigger to generate the pseudo-random number subset for each modal force. 
The means and standard deviations of these pseudo-random subsets are 
specified by the user as AV(1), and SD(1), and AV(2) and SD(2), respec- 
tively. 

This logic assures that the same pseudo-random number subset 
will be used for each case with respect to Card Form 4 parameter varia- 
tion. 

Block C . After the construction of the random sequence component, func- 
tions are mixed. The equations of mixing are given by: 

F(L) = R(1 , 1) * F(L) + R(l,2) * G(L) , 

t - i ? rvi _o 

xu ~ x } L. v/h- -J y 

G(L) = R(2,l) * F(L) + R(2,2) * GfL) , 

where F is the forcing function for the 1st mode and G for the 2nd 
mode. Obviously, no mixing occurs if 

R(l,2) = R(2,l) = 0.0 in the control parameters.. 

The real and imaginary components of the resulting functions 
may be plotted against frequency for each segment number which is an 
integer power of 2. These plots are specified values of ISRS and IERS 
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on control Card Form 3 upon the user's request. The subset of each 
segment from point No. ISRS to IERS will be plotted. 

Block D . A cyclic component is added to the pseudo-random subset con- 
tained in QR1(L) and QR2(L) as follows: 

QR1(L) = QR1(L) + C(l) * COS * NP(1) * NTM + TH(1)] 

for L=1 , 2049 and similarly for QR2(L). NPC1) and NP(2) are the fre- 
quencies in cycles/rad. TH(1) and TH(2) are the associated phase dis- 
placements. These parameters are specified by the user on Card Form 4. 

In addition to adding cyclic components to the forcing func- 
tions, positive gate functions can also be added. The gate function 
generator, called DST in the program, is a function of NTM. The first 
'rise* points of the functions are given by LD(1) for Mode 1, and LD(2) 
for Mode 2, on Card Form 1 . These point No.s refer to the value of NTM„ 
not the segment index. The number of 'rise' points is given by 
UJ(1) and UW(2). The height of each gate is unity. 

Block E . The forcing function and initial conditions are complex -mul- 
tiplied by the modal frequency response functions. 

The real parts are extended symmetrically about CF(2049) or 
CG(2049) respectively. The imaginary parts are also extended but anti- 
symmetrically. IM(CF(2049)) and IM(CF(2049)) are set to zero. The re- 
sulting complex functions, of length 4096 each, are inverse transformed 
into the time domain. The imaginary components of the resulting time 
series are zero-valued. The real valued time series are modal velocity 
data. 
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Block F . Initial conditions for the second and subsequent segments of 
simulation are determined by the procedure below. 

The PKPLT subroutine is called to find the final peak in the 
modal velocity signals of each segment. These peak velocities are the 
respective initial velocities of the following segment. 

The modal displacements are zero when the modal velocities peak 
The initital displacements for the second and all subsequent segments 
are, therefore, zero. 

Block G. STN Is the value of signal to noise ratio which is an option 
to include the random noise to the total system response. The noise 
signals can be obtained by RANDN subroutine if STM is not equal to 
zero. The system response is, then, simply the summation of simulated 
response signals and the noise. The standard deviation of the noise 
component Is equal to the absolute maximum value of the response di- 
vided by STN which is established by input. 

The resulting data segments are input data for the RNDMC1 
algorithm (Randomdec) . The Fortran array SIG(L), L=l, 4096 is the 
sum of the two modal response segments and is used by RNDMC1 to 
establish trigger points, i.e., points which become origins of the 
ensemble summing operation which generates the Signature. The Fortran 
arrays RF (L) , and RG(L), L=l, 4096 are used to hold the components of 
the Signature which arise within the current segment of Modes 1 and 2 
respectively. On completion of RNDMC1 algorithm, RF(L) is summed in- 
to SRF CL) and RG(L) is summed into SRG(L), the array S(L), L=l, 4096, 
is the cumulative sum of the two modal Signatures. YRMS is the RMS 
value of SIG (X) , L=l, 4096. 
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The following diagrams show the logic of the RNDMC1 algorithm 
and the application of its related control parameters, from Card 
Form 2. 





A- 20 









:n 


K=K+1 

I- J 


< N=1 , NPRTy 


RF (N) = RF(N)+F(1,N+1-II) 
RGfN] = RG (N)+F(2 jN+1-11) 
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4.2 Curvft Procedure 


The following diagrams show the Curve-Fit algorithm. - 










CF Procedures 


The execution path of the procedure is shown in the flow diagram 
below. Contents in lettered blocks in the flow diagi'am are further de- 
scribed in the correspondingly lettered paragraphs following the flow diagram. 
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Appendix B 
SAMPLE TEST CASE 


PROGRAM RADCUF 




INPUT DATA 


(Refer to Pages A- 5 to A- 10) 
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SIMULATED RESPONSE SIGNAL C MODE 13 


2fi.i9 



■9 < TaS5I« .0333S *500=103 3 


(Set NSF to 1 to Expand Graph to Occupy Full Page) 
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m/m 



J03 WO 60232 FASE 3 

HQDE 1 SIGNATURE 


ss.s 
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lU/Bfft 


JG3 ua 60232 P9SE * 


MODE 2 SIGNATURE 


i.e 
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IN/SffG 



JOS wo 60232 P«SE S 

SYSTEM SIGNATURE 


gs.s 
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PEAK PLOT OF MODE 1 SIGNATURE FOR 
SIMULATION SIGNAL PARAMETERS 


FCN 

FORM 

FREQ 

RANO0H 

BOXCAR FORCES 

•VJAC 1 

3 a 

21<40 

NUC 1 3 

- 535.0 

PCI. 13 «*0.050 

LOC 1 3 

- 9 . 99999 X 10 

NAC23» 

2190 

NUC 2 3 

« S&S.O 

PU. 25 «0.0!0 

L_ W C 1 3 

a 

NS 

«B 

0 

DAMPING 

PCS, 13 “0.050 

UOC23 

• 9 . 99999 X 10 

NC 

OB 

0 

zcn 

« 0.0317 

PC2, S3 “0.000 

LWC2 3 


INIT 

COND 

ZC23 

« 0.0053 

AVC 1 3 •0.000 

CYCLIC FORCES 

PC 1 3 

as 

0.000 

OTHER 

©DC 1 3 “1.000 

CC 1 3 

“ 0.000 

Qcn 

■3 

98.000 

D OMGA*j 0. 1257 

AVC s 3 «0.000 

NIP 

• 0.100 

P£23 

«a 

0.000 

AC 1 3 

- 29.6 

BDC23 “ 1.000 

TWC 1 3 

- 0.000 

CC23 

83 

C .000 

ACS 3 

— 93-3 


CC23 

• 0.000 




muc a 3 

- 0.0230 


N2P 

•0.100 




HUC23 

- 6.0520 


THC23 

• 0.000 


SEGMENT i 
RAND0MDEC PARAMETERS 

METHODS e 
NO. OF SEGMENTS i* 
TOTAL TIMES TRIGGERED 23* 
NO. OF URNS ON AV6 SIG i 


MAXIMUM PEAK VALUE = 


2.85043310 


•*03 



B IS 28 28 93 58 68 78 S8 58 183 

PEAK NO. 
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PEAK PLOT OP SYSTEM SIGNATURE PGR 


FCN 

FORM 

FREQ 

RANDOM 

BOXCAR FORCES 

MAt 13™ 21*58 

►auc a 3 

« 535.0 

Rtl.I 5 =Q.G5Q 

UOC 1 3 

- 5.95333X10 

NA C 2 3 =» JHfi 

fcaus 2 3 

«= 361.0 

RC1.2J =>0.010 

LW£ 1 3 

n 

M3 

® s 

DAMPING 

RCS. 13 >=0.050 

LOC23 

« 9.999395I0'** 5 

NC 

aa 5 

2U3 

« 0.0317 

0(2. 22=0.000 

Ltft2 3 

m 

IN IT COMO 

ZC 2 3 

« 0.0058 

AVC 1 S =0.000 

CYCLIC FORCES 

ptn 

-» 0.000 

OTHER 

eoc i a = 1.000 

ccn 

—c.ooc 

CC13 

® 5,3.000 

O OH£A 

■» 0.5357 

AVC 2 3 =0.000 

NIP 

-0. 100 

P«2D 

0.000 

AC a 3 

- 2*5. fc 

00(2 3 =1.000 

TWC 1 3 

-0.000 

oca: 

-93.000 

AC 2 3 

— 93.3 


CC23 

-*0.000 



MUJ a 3 

- 0.0240 


*3 25=- 

«®0. ICO 



MiJC 3 3 

«■ 0. 0520 


THC2 3 

—0.000 


SEGMENT i 
RANDGHDEC PARAMETERS 

METHODS 0 
NO, OF SEGMENTS i* 
TOTAL TIMES TRIGGERED 23« 
NO. OF ITRNS ON AVG S15 i 


MAXIMUM PEAK VALUE 


2.9332SI10 


+03 



a is £0 


£0 « SS 60 

PEAK N0. 


70 S3 


100 
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JOS NO 60232 PAGE 61 

INPUT SIGNATURE TO BE CURyE-FITTED 


SS-C 
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JOB m 6G232 PASS *3 

COMPUTED VALUE GF LQG t HE AN-SQUAREO ERROR 3 FOR EACH ITERATION 
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PJcS <W 


CONVERGENCE PATH OF A( 13 
12 .*! 30 



J?33vnBii 08. 










IN/SffC 


JOB KO 60232 


PflSE 


CONVERGENCE PATH OF AC 43 


17 


-. 3 $$ 
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JC3 W3 60232 P»S£ <!3 

INPUT SIGNATURE TO BE CURVE-FITTED 


za.e 
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iti/atrc 


JOB MO 60232 


PASS 99 


LEAST-SQUARE APPROXIMATION ON ITERATION 13 


».e 
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JOS m fc0232 F£5E 50 

COMPUTED VALUE OF LOG t MEAN-SQUARED ERROR 3 FOR EACH ITERATION 


*.e 
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JOB MO 40232 PASE 51 

CONVERGENCE PATH OF AC U ' 


11. 6M 



I'5at3a73e^ no. 
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JOB »0 60S32 


P3SE 53 


CONVERGENCE PATH Qf AC 2) 
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JOB HO 6023B PASS 53 

CONVERGENCE PATH OF AC 33 


.89780 
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JOB «IO M232 'P? SE 5*4 

CONVERGENCE PATH OF AC 43 





3 <5 5 4 7 S f 10 11 12 

ITJSatlOa as. 
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JOB 30 


60232 


PAS* 55 


CONVERGENCE PATH GF AC 53 
2.0 
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940/&SC 


J03 m 60232 MS* 56 

CONVERGENCE PATH OF AC 63 



ITBMAVICP MS. 


BA0/8&C 


J03 SO 


60332 


PISE 57 


CONVERGENCE PATH QF AC 73 
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